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Anaerobic digestion enables the water industry to treat wastewater as a resource for generating 
energy and recovering valuable by-products. The complexity of the anaerobic digestion process 
has motivated the development of complex models. However, this complexity makes it intractable 
to pin-point stability and emergent behaviour. Here, the widely used Anaerobic Digestion Model 
No. 1 (ADM1) has been reduced to its very backbone, a syntrophic two-tiered microbial l food 
chain’ and a slightly more complex three-tiered microbial ‘ food web’, with their stability analysed 
as function of the inflowing substrate concentration and dilution rate. Parameterised for phenol and 
chlorophenol degradation, steady-states were always stable and non-oscillatory. Low input concen¬ 
trations of chlorophenol were sufficient to maintain chlorophenol- and phenol-degrading populations 
but resulted in poor conversion and a hydrogen flux that was too low to sustain hydrogenotrophic 
methanogens. The addition of hydrogen and phenol boosted the populations of all three organisms, 
resulting in the counterintuitive phenomena that (i) the phenol degraders were stimulated by adding 
hydrogen, even though hydrogen inhibits phenol degradation, and (ii) the dechlorinators indirectly 
benefitted from measures that stimulated their hydrogenotrophic competitors; both phenomena hint 
at emergent behaviour. 


I. INTRODUCTION 

Microbial degradation of organic compounds in 
methanogenic environments is a sequential process catal¬ 
ysed by a series of different micro-organisms. Syntrophy 
plays a pivotal role in these feeding webs: degradation 
of compounds like propionate and phenol is only sustain¬ 
able if their degradation products, hydrogen and acetate, 
are removed by methanogens. The thermodynamic ra¬ 
tional behind syntrophy is well understood, but its ki¬ 
netic framework is less established. This raises questions 
about the stability of these feeding chains and the fac¬ 
tors that govern them. As a first step towards answering 
these questions a simple mathematical model was pre¬ 
viously developed describing the interactions in a two- 
tiered feeding chain, populated with a set of parameters 
that apply to propionate degraders and hydrogenotrophic 
methanogens [T; . Mathematical analysis of the model in¬ 
dicated that the system was always stable: there were no 
conditions where the populations of the two organisms 
oscillate or show other forms of emergent behaviour. 

The objective of the present paper is to introduce an 
additional organism into a similar feeding chain and eval¬ 
uate its effect on stability and potential emergent be- 
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haviour of the resulting ‘food web’. The organism of 
choice is a chlorophenol-dechlorinating bacterium. The 
other two organisms are a phenol degrader and a hy¬ 
drogenotrophic methanogen. The complete removal of 
phenolic compounds from the system is hereby referred 
to as chlorophenol mineralisation. The salient feature of 
the chlorophenol degrader here is that production of phe¬ 
nol is coupled to consumption of hydrogen by hydrogen 
cycling. Thus, as a hydrogen consumer, the dechlorina- 
tor competes with the methanogen for hydrogen m®. 
The working hypotheses are (i) that the dechlorinator 
can (partially) replace the methanogen as the syntroph 
in a phenol-degrading consortium, and (ii) that introduc¬ 
tion of this organism can potentially lead to unexpected 
emergent behaviour related to the intricacies of the multi¬ 
species relationships. 

It has been shown that deterministic modelling of bio¬ 
logical systems, typically through a system of coupled or¬ 
dinary differential equations (ODEs), provides important 
understanding of these often complex processes, specif¬ 
ically in determining changes to the system behaviour 
given perturbations in the inputs. For anaerobic diges¬ 
tion, higher dimensional models are useful for capturing 
the phenomenological behaviour of the multi-step pro¬ 
cesses and are often the de-facto method for understand¬ 
ing plant operation [4|. On the other hand, simplified or 
reduced models have received more attention in process 
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monitoring, control design 0 0, and optimisation [7]. 
Simplified models have also been applied to determine 
the global 00 or local m stability of the system un¬ 
der investigation. Typically, these are two-species mod¬ 
els, whereas here we present a case considering a three- 
species food web. 

Whilst a strictly analytical approach is not possible 
given the dimensions of the ensuing model, the approach 
taken here derives analytic expressions for all steady- 
states, supported by numerical simulations to determine 
the regions of local stability within sensible operating 
conditions. As such, it is possible to gain a formal un¬ 
derstanding of the emergent properties of all states. 
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II. MODEL AND METHOD 


Here the functions /o, /i and fa are defined as 


A. Mechanistic model of a three-tiered food web 


The model developed here is based on Anaerobic Diges¬ 
tion Model No. 1 (ADM1) 0. The general model here 
has six components, three substrate and three biomass 
variables, from which a sub-model describing phenol 
degradation and the extension of the full model to in¬ 
clude addition of extraneous substrates, are formed and 
described in the relevant sections. The chlorophenol 
degrader utilises both chlorophenol and hydrogen for 
growth, producing phenol as a product. Phenol is con¬ 
sumed by the phenol degrader forming hydrogen, which 
also is inhibitory to its growth. The methanogen scav¬ 
enges this hydrogen and acts as the primary syntroph. 

The time-dependent substrate and biomass concentra¬ 
tions are denoted S(t) and A(t), respectively, with sub¬ 
scripts c h, P h and h 2 referring to chlorophenol, phenol and 
hydrogen components, respectively. The growth func¬ 
tions are of Monod form with the inclusion of a product 
inhibition term, A/ h 2 - Chlorophenol, phenol and hydro¬ 
gen are introduced with an input concentration .S/h.in, 
'Sph.in, *^H 2 ,in 5 respectively, und U dilution rute D. The 
inhibition of hydrogen on the phenol degrader is defined 
as 
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where Ksm 2 ,c and A's.h 2 are the half-saturation con¬ 
stants for hydrogen in the chlorophenol degrader and 
hydrogenotrophic methanogen, respectively. The values 
used for the various parameters are listed in Table [T] and 
their derivation is presented |A| The value 224/208 rep¬ 
resents the fraction of chlorophenol chemical oxygen de¬ 
mand (COD) converted to phenol, 32/224 is the fraction 
of phenol converted to hydrogen, and 16/208 is the frac¬ 
tion of hydrogen COD consumed by A' c h- 


B. Dimensionless form 

It is beneficial to scale unit-dependent equations to a 
dimensionless form; this significantly reduces the number 
of parameters describing the dynamics, thereby simplify¬ 
ing the subsequent analyses. Using the notation of Xu 
et al. [I], after Baltzis and Fredrickson m. the following 
dimensionless terms are defined 
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The substrate and biomass concentrations then evolve 
according to the six-dimensional dynamical system of 


ODEs 
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With these transformations, the dynamical system of 
ODEs given in Equations reduces to 
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Parameters Nominal values Units 
km,ch 29 

K s , ch 0.053 

Uch 0.019 


fcm, ph 

26 

Ks, ph 

0.302 

T ph 

0.04 

^m, H2 

35 

A’s,h 2 

2 . 5 x 10 " 

A's,H 2 ,c 

1.0x10"' 

Th 2 

0.06 

^dec,i 

0.02 

Ki,h 2 

3 . 5 x 10 "' 


kgCODs/kgCODx/d 

kgCOD/m 3 

kgCODx/kgCODs 

kgCODs/kgCODx/d 

kgCOD/m 3 

kgCODx/kgCODs 

kgCODs/kgCODx/d 

kgCOD/m 3 

kgCOD/m 3 

kgCODx/kgCODs 

d" 1 

kgCOD/m 3 


TABLE I. Parameters used in the two-tiered food chain and 
three-tiered food web models, where d represents days, and 
CODx and CODs are the Chemical Oxygen Demand (COD) 
of the biomass and substrate, respectively. The derivation of 
these parameter values is given in [A] 
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C. Steady-states 


The steady-states of this system of ODEs are obtained 
by setting g, = 0 (for i = 1,..., 6) and solving simulta¬ 
neously. This predicts eight possible steady-states: 

SSI:: The trivial steady-state where all three popula¬ 
tions are washed out (xo = x\ = X2 = 0). 


SS2:: Only the methanogen population is maintained 

(x 0 = x 1 = 0 ,x 2 ^ 0 ). 
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SS3:: The phenol degraders and methanogens are 
washed out (xo 7^ 0, X\ = X 2 = 0). 

SS4:: The hydrogenotrophic methanogens are washed 
out while the chlorophenol and phenol degraders 
are maintained (xq 7^ 0, X\ 7^ 0, X 2 = 0). 

SS5:: Only the phenol degraders are washed out (xo 7^ 
0,xi = 0, x 2 7^ 0). 

SS6:: All three populations are present (xo 7^ 0,Xi 7^ 
0, x 2 7^ 0). 

SS7:: Only the phenol degraders are present (x 0 = 
0, xi 7^ 0, x 2 =0). 

SS8:: The phenol degraders and methanogens are main¬ 
tained in the system (xo = 0, Xi 7^ 0, x 2 7^ 0). 

For full chlorophenol mineralisation and, as such, the 

only desired operating condition, SS6 must be stable. 


D. Stability analysis 

where, for simplicity, the following parameters are intro- The challenge of analytically characterising dynami- 
duced cal systems of ODEs is well understood. Typically, for 
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systems of high dimensions, it is intractable to finding 
explicit solutions and one must resort to obtaining so¬ 
lutions numerically. For example, the Routh-Hurwitz 
theorem allows for an explicit analysis of the stability 
of steady-states, but is intractable beyond five dimen¬ 
sions [ 12 ] and is excluded here. Relying on numerical 
solutions, however, is problematic as the extent and res¬ 
olution of the results are limited by the choice of model 
parameters and the computational resources. Neverthe¬ 
less, here analytical expressions may be found for some 
of the steady-states that can inform general rules about 
their viability. To present the results of the stability anal¬ 
ysis coherently, all steady-states are solved numerically 
for a realistic range of operational and kinetic parameters 
with a suitable resolution of 5000 x 5000 solution points. 
A brief discussion of the numerical methods used in this 
analysis is provided in|B| 


In order that a steady-state be meaningful, all vari¬ 
able concentrations must be non-negative, while Eq. i> 
also gives the condition that S c h < Sch.iiu ( or so < it/ in 
dimensionless form). For model extensions with phenol 
and hydrogen addition, the conditions that S p h < S p h,in 
(si < u g ) and Sh 2 < Sh 2 ,in (s2 < Uh) are also necessary. 
It is well established that the stability of a system of au¬ 
tonomous ordinary differential equations (ODEs) can be 
determined by investigating the eigenvalues of the cor¬ 
responding Jacobian matrix m ■ The Jacobian for the 
system of Equations ( TT]JT6 ) corresponds to the 6x6 ma¬ 
trix 



FIG. 1 . Steady-state diagram for operational parameters D 
and Sph.in, in the two-tier phenol model. Green lines be¬ 
tween steady-states: explicit functions of the steady-state 
partitions, after Xu et al. [TJ. 


III. RESULTS 
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Note that for the two-tier model, for which X c h = 
S c h = 0 , the Jacobian reduces to a 4 x 4 matrix. The 
Jacobian is then evaluated at a given steady-state, de¬ 
noted Jssi (i=l,...,8), and its eigenvalues calculated. If 
the real parts of all the eigenvalues obtained from the 
Jacobian are negative, then the state is stable. If one or 
more of the eigenvalues have a positive real part then the 
steady-state is unstable. 

In Section m a thorough investigation is carried out 
in order to determine when each steady-state is viable. 
Firstly, a two-tier model describing only phenol degra¬ 
dation is considered, where Xq = sq = 0, reducing the 
system of ODEs to four equations. Subsequently, the 
full model for chlorophenol mineralisation is studied. 


A. Two-tier phenol model 


An initial analysis was performed on the two-species 
model feeding on phenol. As expected from the model 
structure, the results are similar to those shown for a 
propionate-degrading bi-culture with a maintenance term 
(kdec) [I], where three stable steady-states emerge: SSI, 
SS 7 and SS8. Figure [I] shows the steady-state diagram 
for the phenol model, demonstrating that the system has 
three mutually exclusive stable steady-states under the 
range of parameters and operating conditions chosen. 
Numerical simulations under conditions producing SS 7 
(not shown) indicated that the concentration of phenol 
degraders are comparatively low (compared to SS8 pop¬ 
ulations), resulting in low hydrogen production. This 
in turn results in the washout of the methanogen, whilst 
the phenol degraders, under significant phenol concentra¬ 
tions, can be maintained up to dilution rates equivalent 
to their theoretical maximum growth rate minus the de¬ 
cay constant, /i m ax,ph — ^dec.ph = 1 . 02 d _1 . Although nu¬ 
merical analysis of the steady-states was performed here, 
it is possible to get explicit quadratic functions of the pa¬ 
rameter pair (S'ph/n, D) for the steady-state partitions. 
Using the method described in Xu et al. T], the functions 
for the partitions are (in non-dimensionless form) 
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where Fi = 0 relates to methanogen washout and 
F 2 = 0 to phenol degrader washout. These partitions 
are shown by the green lines in Fig. [T] 

Having observed that the phenol model is always sta¬ 
ble and the desired operational state (SS8) is bounded 
asymptotically at higher phenol input concentrations, the 
next sections explore the properties of the model when 
extended to a three species dechlormating tri-culture 
with hydrogen cycling. 


B. Three-tier chlorophenol model 


Solutions were found for the three steady-states iden¬ 
tified (SSI, SS 4 and SS6). Defining 


E = 
H = 


duo 

ds 0 

dm 

ds 2 


> 0 ; F = 


< 0; I = 


d^o 

ds 2 

d ^ 2 

ds 2 


> 0 ; 

> 0 , 


dni 

ds 


( 20 ) 


then the Jacobians for SSI, SS 4 and SS6, respectively, 
are 
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A full analysis of the stability for each steady-state 
is provided in[C] By inspecting the eigenvalues (A,;, i = 
1 ,..., 6) of the Jacobian matrix for each steady-state, the 
following conditions for stability are found: 


SSI:: Meaningful and stable, with so = Uf. 


SS4:: Meaningful when the roots of the cubic function 
(Eq. ( |C 11 [ )) for s 0 are positive and the conditions 
a < (so/(l+so))—and a < ( 0 i/(l +K I s 2 ))- 
ks are met. Additionally, it can be seen from 


Eq. (C 18 ) that Xq can only be meaningful if wqoji > 


uj 2 . Solving this inequality results in the condition 
( 1 —lch)(l—Fph) > 0 . 5 , such that Xq is meaning¬ 
ful provided 0 < Yph < 0.5 and 0 < F c h < 0 . 5 . 
Stability requires that a > —kc and R(A,;) < 0 
(i = 2.6). 


SS6:: Meaningful for the following conditions 
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and stable when R(Aj) <0 (i = 1 ,..., 6). 


Numerical simulation of the model was performed us¬ 
ing the parameters listed in Table [TJ over a range of oper¬ 
ating conditions ( a,u /). The stable steady-state regions 
are shown in Fig. [ 2 ] As determined from the stability 
analysis, SSI is meaningful under all operating conditions 
and bistable steady-states with SS 4 and SS6 are found in 
specific regions within the numerical limits of the simula¬ 
tion, depending on the initial conditions. An interesting 
phenomenon not observed in previous two-tiered models 
is the existence of SSI at low dilution rates but increasing 
chlorophenol input. It is likely that the condition Uf > Sq 
is not satisfied in this region and suggests that the system 
does not produce enough metabolites to sustain the syn- 
trophic populations at such low substrate concentrations 
and flow rates. 


C. Three-tier chlorophenol model with hydrogen 
addition 

Given that the model under analysis includes the pos¬ 
sibility of adding extraneous substrates (phenol and hy¬ 
drogen) to the system, supplementing the flux produced 
by the biomass, it is fairly straightforward to extend its 
mathematical analysis with these terms included. It is 
hypothesised that the addition of hydrogen to the sys¬ 
tem will result in the retention of the methanogen at 
dilution rates up to its maximum growth rate (/r max ). 
It is anticipated that this leads to the extension of full 
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S ch,in( k 9 COD/m ) 


FIG. 2. Steady-state diagram for operational parameters 
D and S' c h,in in the three-tier chlorophenol model (.Sph.m = 
S , H 2 ,in = 0). 


chlorophenol mineralisation by allowing the establish¬ 
ment of a methanogenic population that can be main¬ 
tained given the availability of hydrogen in high enough 
concentrations. 

When hydrogen addition is included, the hydrogen in¬ 
put term, 1 S , H 2) i n > 0 in Eq. Q, and the dimensionless 
form, Uh > 0 in Eq. (16). A further three steady-states 


are defined with this model; SS2, SS3 and SS5. 

Following the same approach as for the three-tier 
chlorophenol model, six Jacobian matrices were derived 
for each steady-state, and their eigenvalues A* (where 
i = 1,..., 6) found. In this case, the conditions for sta¬ 
bility for each steady-state are as follows 


SSI:: Meaningful and stable providing 


Uf u h 
1 +Uf Kp+Uh 


—Ua < a; 


<^2 Uh 

1 +u h 


—kc < a. 


Given that, in this model, kA = kc, then for 
(f >2 > 1, the second condition becomes the pre¬ 
vailing factor in determining stability, and can be 
reduced to Mmax.H 2 ~kc < a, where Mmax, h 2 = 

^H 2 /^m,ch^ ch- 

SS2:: Meaningful when a < </> 2 — fcc and Uh > 

S 2 and stable when K p{u f+l){cj )2 — kc)+kc~Uf+ 
kcuf/(Kp—l)(uf+l) < a. Here it can be noted 
that the stability of SS2 does not depend on Uh , 
and as such will remain fixed for any hydrogen in¬ 
put concentration. 


SS3:: Meaningful given the roots of the quadratic func¬ 
tion (Eq. (C3)) for So are positive, Uf > So and 
Uh > ^(uf—So). For steady-state stability, it is 
necessary that a > Mi —kp, ol > M 2 —fcc and for 
R(Ai) < 0 [i = 4,5,6). Here, the stability of SS3 


is dependent on both it/ and Uh■ The first condi¬ 
tion dictates that at higher it/, a greater amount 
of Uh is also required to maintain stability, whereas 
the second condition indicates that at a fixed it/, 
increasing Uh will increase the dilution rate under 
which SS3 is stable up to a boundary described by 
a function /(«/, Mmax,H 2 -fcc)- 


SS4:: Meaningful given the roots of the cubic function 
(Eq. ( |C11[ )) for so are positive and the conditions 
a < (so/(l+so))-k A , a < (</>i /(l+K I s 2 ))-k B , 
Uh > W 1 S 1 +S 2 and x\ > 0 are met. Stability 
is assured when a > M 2 ~kc and for R(A i) < 0 
(i = 2,..., 6). As Uh affects So through the highly 
non-linear cubic expression, it is difficult to draw 
any conclusion about its effect on stability by ana¬ 
lytical means. 


SS5:: Meaningful when the conditions <f >2 > a+kc, 
s 2 > (K P +s 2 )(a+kA), Uf > s 0 and u h > s 2 + 
ui 2 +{uj 2 kAXo/a) are met. Stability of this steady- 
state is assured when a > fii—kp and R(Ai) < 0 
(i = 3 ,..., 6). 

SS6:: Meaningful when a < (j> 2 —kc, ct < 

(4>i/{l+KiS 2 ))-k B , a < {4>i/(K P +S2))-k A , 

Uf > So, a < io^kAXo/{ si—ujqXq) and a; 2 > 0. 
SS6 is stable when R(Ai) < 0 (i = 1,..., 6). The 
analysis shows that SS6 stability is affected by 
hydrogen addition. 


The results from the stability analysis are visualised in 
Fig- i for different concentrations of hydrogen addition. 
As can be seen, four steady-states (SS3-6) are affected by 
increasing concentrations of hydrogen. At very low hy¬ 
drogen concentration addition (Fig. [ 3 ] (a)), SS3 is stable 
in the region between (0,0) and SS2 (SSI is not shown as 
it remains fixed at D = M maXi H 2 —fcd ec = 2.08d~ 1 , where 
A/max is the non-dimensionless form of Mmax)- Moreover, 
at lower dilution rates, the system has two observable 
bistable conditions: 1) between SS3 and SS6 at the low¬ 
est dilution rates; and 2) between SS3 and SS4, before 
SS3 becomes the single stable steady-state at higher di¬ 
lution rates. This can be rationalised as follows: 


• The chlorophenol degrader can outcompete the 
methanogen for the additional hydrogen in the sys¬ 
tem. Depending on the initial biomass concentra¬ 
tions, this can result in complete washout of the 
phenol degrader and methanogen at low dilution 
rates. However, at higher chlorophenol concentra¬ 
tions, the system can be in stable equilibrium with 
all organisms present; 

• As the dilution rate increases, the methanogen is 
washed from the system as the chlorophenol de¬ 
grader outcompetes it to become the sole hydrogen 
scavenger, as seen in Fig. [2] 
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• Eventually at higher dilution rates, the phenol de¬ 
grader is washed out and full chlorophenol mineral¬ 
isation is no longer possible. Nevertheless, beyond 
the theoretical M max c h, the methanogen is able to 
grow again (SS2), as there is a constant supply of 
substrate without an active competitor. 

With increasing hydrogen addition interesting phe¬ 
nomena are observed (Fig. [ 3 ] (b)-(d)). Firstly, the bista¬ 
bility between SS3 and SS6 is replaced by a single steady- 
state, in which all organisms are present. Similarly, the 
SS3-SS4 bistability is replaced by only SS4, although 
this requires greater hydrogen concentrations at higher 
chlorophenol input. A new steady-state also begins to 
emerge (SS5), starting at low chlorophenol input, which 
eventually subsumes the other two varying steady-states 
(SS3 and SS4) at high hydrogen concentrations. The un¬ 
derstanding from this is that at dilution rates higher than 
the upper boundary of SS6 and below the lower bound¬ 


ary of SS2, the phenol degrader is naturally washed out, 
but the presence of increasing amounts of hydrogen al¬ 
lows the methanogen to exist whilst in competition with 
the chlorophenol degrader. 

Most importantly, the extent of SS6 is also observed to 
increase below a specific <S' c h,in, with the additional hy¬ 
drogen resource stabilising the methanogen population 
at higher dilution rates, limited now only by its maxi¬ 
mum growth rate (M maXi H 2 ). Above this dilution rate, 
as with the three-tiered model, the production of phenol 
is limited by the stoichiometry of chlorophenol degra¬ 
dation. This restriction results in a limited supply of 
phenol under these operating conditions, which cannot 
produce enough growth in the phenol degrader to sus¬ 
tain its population. However, the increasing abundance 
of additional hydrogen leads to the establishment of a 
two-species trophic level with no intermediary organism 
(SS5). Nevertheless, the ability to extend the stable re- 




















gion of SS6 with addition of hydrogen is of note, partic- D. Three-tier chlorophenol model with phenol 

ularly for practitioners. addition 


It has been shown that the inclusion of hydrogen addi¬ 
tion in the three-tiered system can lead to an increase in 
the stable region for full chlorophenol mineralisation with 
regard to standard operating parameters. However, this 
is limited to a defined operating range. At lower S' c h ) i n 
there is washout of phenol degraders and establishment 
of methanogens at increasing hydrogen concentrations. 
It appears that a form of competitive exclusion occurs at 
certain (D,S C h,i n )- However, unlike the classical case, 
here the exclusion principle dictates that consumer A 
(A ch ) is never excluded and either consumer B (A H2 ) 
or producer A (X p h) are washed from the system. This 
can be seen in the system simulation shown in Fig. II 
which considers a single operational point (D = 0.36d , 

S c h,in = 0.30kgCOD/m 3 ), with increasing Sn^m- 


For such a non-linear system the addition of hydrogen 
could lead to deleterious and potentially unexpected be¬ 
haviour, especially given its inhibitory effect on the phe¬ 
nol degraders. It is therefore useful to examine the effect 
of adding phenol to the system, both as a sole exoge¬ 
nous substrate, but also in combination with hydrogen 
addition. In this way, it may be determined under what 
conditions, if any, full chlorophenol mineralisation can be 
extended in relation to the two operational parameters 
(D and S'ch.in)- 

When phenol is added to the system the phenol inflow 
term, SV,h,in from Eq. ([6]) and u g from the dimensionless 
Eq. (15), become positive. An additional two steady- 
states are defined with this model, SS7 and SS8. 

The conditions for stability in each steady-state, deter¬ 
mined in the same manner as previously described, are 
as follows: 


At low hydrogen addition, the concentration is 
not enough to sustain the methanogen population 
whilst a small population of dechlorinators are present 
(6x 10~ 5 kgCOD/m 3 ). The hydrogen not utilised by the 
chlorophenol degrader reaches an equilibrium that in¬ 
hibits the phenol degraders and they are washed out of 
the system and, thus, SS3 is the dominant steady-state 
under these conditions for the initial conditions tested. 
As hydrogen addition is increased, an interesting prop¬ 
erty of the system emerges. Again, the chlorophenol de¬ 
grader is able to utilise the extra hydrogen to produce 
more biomass and more phenol. However, the concentra¬ 
tions of additional hydrogen are not enough to sustain 
the methanogen population and the equilibrium concen¬ 
tration returns to a value approximate to the previous 
case. Whilst this exerts an inhibitory effect on the phenol 
degraders, the additional phenol substrate availability is 
enough to maintain a population in the system (steady- 
state not shown in Fig. [ 4 ]) , which switches the system to 
SS4. In the final case, hydrogen addition is further in¬ 
creased, allowing chlorophenol degraders to utilise most 
of the available chlorophenol and reach close to its maxi¬ 
mum growth. As this growth is asymptotic, and because 
of the surplus hydrogen available, the methanogens are 
now able to utilise the remaining concentrations to main¬ 
tain their population in the system. This results in an 
increased hydrogen equilibrium that exerts a greater in¬ 
hibition on the phenol degraders, to such an extent that 
the additional phenol availability is not enough to avoid 
washout and the system moves to SS5. 


It can be seen that this situation is fairly complex and 
is dependent on a number of factors. Principally, the 
equilibrium concentration of hydrogen is determined by 
the dilution rate, which in turn limits the growth of either 
producer A or consumer B. 


SSI:: Stable for c/)iu g /(l+u g )—k p < ot. For high con¬ 
centrations of phenol, this is equivalent to /^ m ax,ph— 

ks < Oi, where /Tmax.ph = ^m,ph1 pli /^"m.ch 7 c:h - 

However, at lower phenol concentrations, the effect 
of phenol lowers the minimum dilution rate under 
which SSI is stable. 


SS2:: Never stable as a contradiction occurs in the 
terms describing S 2 and x 2 . For these variables 
to be meaningful, the conditions a—kc < (f> 2 and 
a—kc > 4 >2 must be satisfied, which is invalid (see 
Eqs. (Cl and (fC2])). 


SS3:: Never stable as a contradiction occurs in the terms 
describing s 2 and Xq. With Uf, =0, So > it/ for 
S 2 to be positive. This condition, however, results 
in Xn always being negative and meaningless (see 
Eqs. (Cl) and (Si). 


SS4:: Meaningful given the roots of the cubic function 
(Eq. (Clip ) for so are positive and the conditions 
a < (so/(l+so))-k A , a < (^i/(1+A"/s 2 ))-A; b , 
Xq > 0 and x\ > 0 are met. Stability is assured 
when a > ^2 ~kc a nd th e roots of the characteris¬ 
tic polynomial (Eq. (C31)) have negative real parts. 


In the special case that hydrogen is also added to 
the system, it should be noted that if u g < si, then 
the condition Uh > S 2 must be strictly observed. 
Similarly, if Uh < S 2 , then u g must be greater than 
Sl- 


SS5:: Never stable as without Uh then X 2 can never be 
positive (see Eq. (C29)). 


SS6:: Meaningful when a < <p 2 ~kc, a < 

(</>i/(1+A7s 2 ))-£:b, a < (4>i/(K P +s 2 ))-kA, 

iif > so, xq > 0 and a < (uiiXi— uqa^o)/^- With 
these conditions met, stability of SS6 is guaranteed 
when the roots of the characteristic polynomial 
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FIG. 4. Evolution of three-tier system with H2 addition, where D = 0.36 d _1 , S c h,in = 0.30 and Sj^.in = 2.67xl0 -4 (green), 
5h 2 ,in = lxlCP 1 (red), Sb^m = 2.67 (blue). Initial conditions: X,-i t = 10 —3 , A' p h = 6.5xl0 -5 , Ah 2 = 10 —3 , S c h = 10 -2 , 
S'ph = 5xl0 -2 , Sh 2 = 10 -3 (all kgCOD/m 3 ). 


(Eq. ( |C42 )) have negative real parts. In this case, 
it can be shown that SS6 is influenced by u g . 

SS7:: Meaningful given the ro ots of the quadratic func¬ 
tion for si (Eq. ( |C50 )) are positive and, sub¬ 
sequently, (j)\Si/(l+S\)Ki{a+kB) > 1/Ki, and 
u g > s-|. Therefore, stability is assured when 
a > [iQ-ltA, cr > fj, 2 ~kc a nd th e roots of the char¬ 
acteristic polynomial (Eq. (C53)) have negative real 
parts. It can be seen that the stability of SS7 relies 
on u gi but is also influenced by hydrogen if added 
to the system. 

SS8:: Stable when the following conditions are satis¬ 
fied. Firstly, S 2 is meaningful given (f >2 > a—k c . 
Stability is assured when S 2 is meaningful, a > 
Ho—kA, and the r oots of the characteristic poly¬ 
nomial (Eq. ( C62| ) have negative real parts. In 
general, for SS8 to be stable, a must be smaller 
than /i max ,H 2 and, therefore, Uf should be small. 


Again, the operational parameter plots are used to vi¬ 
sualise the outputs from the stability analysis, as shown 
in Figs [5] (a)-(c). At low phenol addition, emergent 
steady-states (SS6-SS8) are observed at low dilution 
rates, with the bistability present in the standard three¬ 
tiered model comprising the rest of the system below 


the boundary for SSI. As phenol is increased, the two 
bistable regions disappear to be replaced by a single sta¬ 
ble steady-state (SS4 and SS6). Both SS7 and SS8 also 
appear across a greater range of dilutions, although SS8 
is confined to low chlorophenol addition, whilst SS8 even¬ 
tually covers a region contained by the upper bounds of 
SS7 (low S'ch.in) and SS4 (high SA^in), and the boundary 
of SSI. Of greater interest is the observation that the ex¬ 
tent of the desired condition, SS6, is greater for S' c h,i n be¬ 
tween approximately 0 and lkgCOD/m 3 , replacing both 
SS4 and part of SS7. 


E. Three-tier chlorophenol model with bi-substrate 
addition 

Figure [5 (d) shows an example when both hydrogen 
and phenol are added to the system. Here it can be 
seen that the two additional substrates each contribute 
to distinct emergent properties. The hydrogen addi¬ 
tion allows for the maintenance of a methanogen popula¬ 
tion at higher dilution rates, whilst the phenol addition 
maintains the phenol degraders at lower dilution rates. 
Both additional substrates provide the conditions neces¬ 
sary for full chlorophenol mineralisation within a specific 
operational parameter space previously associated with 
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methanogen washout. 

An extension of this analysis was undertaken by ex¬ 
amining the effect of the two supplementary substrate 
additions on the maintenance of SS6 at higher dilution 
rates than observed in the standard three-tiered model. 
Figure [6] shows the steady-state plot for varying concen¬ 
trations of hydrogen and phenol addition for two distinct 
operating conditions observed previously to give rise to 
SS6 with addition of an extraneous substrate. Figure [6] 
(a) demonstrates that both hydrogen and phenol addi¬ 
tion can transform the system from undesirable bistabil¬ 
ity to a stable process with all three organisms present. 
Indeed, it is shown here that addition of either substrate 
can lead to SS6 under relatively high chlorophenol input 
conditions. 

An alternative situation under low chlorophenol ad¬ 
dition is shown in Figure [6] (b). Here, it is clear that 
phenol addition at concentrations above 0.65kgCOD/m 3 


provides the necessary conditions for full chlorophenol 
mineralisation, irregardless of hydrogen addition. Al¬ 
though it can also be seen that hydrogen addition does 
have some effect within a very small range of phenol con¬ 
centration between about 0.63 and 0.65kgCOD/m 3 , al¬ 
lowing for stability of the methanogenic population and 
moving the steady-state from SS4 to SS6, it is relatively 
negligible given the influence of phenol with these oper¬ 
ating parameters. 

IV. BISTABILITY ANALYSIS 

Throughout Section m numerous regions of bistabil¬ 
ity were observed between SS3-SS4 (Fig. [ 3 ] (a), (b) and 
(c)), SS3-SS6 (Fig. H (a) and (b)) and SS1-SS6 (Fig. [TJ). 
In this Section, these bistable regions are investigated in 
order to provide a more in-depth description and under- 
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standing. Through analysis of these results, the steady- 
state having the larger basin of attraction is determined, 
i.e. the most likely outcome for a range of initial condi¬ 
tions. Note that the other four steady-states, even when 
meaningless or unstable, are still present mathematically 
and can impact on the overall dynamics of the system. 

Firstly, the region of bistability between SS3 and SS4 
was considered. Figure [7] shows the final stable steady- 
states (SS3 or SS4) achieved for a range of initial con¬ 
ditions. Note that for SS3 or SS4, = 0 and does 
not deviate from this value throughout, which simplifies 
the analysis required here. When X c h is small, the sys¬ 
tem is likely to head towards SS3 where washout of X p h 
occurs. Otherwise, the system achieves SS4. Very simi¬ 
lar results were obtained for biologically meaningful sub¬ 
strate initial conditions between 1.3xl0 -2 and 1.5xl0 _1 
for S c h, 9.0xl0 -2 and 4.9xl0 _1 for S'ph and 2.2xl0 -12 
and 4.7xl0~ 12 for Sh 2 (kgCOD/m 3 ). These results sug¬ 
gest that as long as sufficient dechlorinating biomass is 



FIG. 7. Phase plot showing final steady-state with initial 
conditions Vh 2 = 0, S c h = 0.14, S'ph = 0.11 and Sh 2 = 5.5x 
10 -5 and parameters D = 0.375d _1 , S c h,i n = 3 and SH 2 ,m = 
0.0267 (All Xi and Si in kgCOD/m 3 ). 


present in the system, the most likely outcome is for both 
Vch and X p h to be present over time. Thus, the basin of 
attraction for SS4 is larger than that for SS3. In addition, 
changing the dilution rate D , S c h,in, or Sh^hi within the 
bistable region once again leads to similar results. 

A similar analysis of the bistable regions for SS3-SS6 
and SS1-SS6 give a more dramatic outcome. If any of the 
three organisms, V p h, V c h or Ah 2 , are initially zero, SS6 
is not obtained. If all of these are present, even in tiny 
quantities, SS6 will be reached. It is known that SSI is 
always meaningful, thus a steady-state always occurs at 
the origin, be it stable or unstable. In the case of bista¬ 
bility, SS3-SS4 or SS3-SS6, the unstable steady-state is 
closer to SS3 rather than SS4 or SS6, respectively. This 
instability thus forces the dynamics away from the total 
washout and towards the state with the least washout. 
This analysis highlights how not only the main parame¬ 
ters of the system can impact the final stable state, but 
also the importance of initial conditions. 


V. DISCUSSION 

In this work, a mathematical analysis of a three-tiered 
‘food web’ comprising three organisms with hydrogen 
addition and inhibition has been presented. Although 
mechanistic models of microbial interactions are some¬ 
what ubiquitous, analytical approaches have been limited 
to simple two-species systems mmm- More recently, 
attempts at analysing more complete systems, such as the 
anaerobic digestion process, have been performed with 
some degree of rigour psum, however the practicality 
of these approaches is limited by the need for generalisa¬ 
tion and increasing numbers of assumptions to allow for 
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mathematical tractability. 

Although the three-tiered system is considered only 
a sub-process of anaerobic digestion and the analysis 
limited to a specific compound, the approach is general 
enough to provide information about the characteristics 
of such a process and allows for the possibility of ex¬ 
tending the work to other systems and other compounds. 
Indeed the results shown here point to the fact that bio¬ 
logical knowledge can inform mathematical approaches, 
whilst mathematical models can indicate or confirm bio¬ 
logical properties of a system, succinctly and rationally. 

Here, four models were analysed by a combination of 
analytical and numerical techniques to obtain informa¬ 
tion regarding the extent and characteristics of system 
stability in a three-organism process anaerobically de¬ 
grading chlorophenol. In this case, only the hydrogen 
part of the pathway was considered, with acetate be¬ 
ing excluded from the analysis. The simple two-tiered 
food chain considering only phenol degradation with two 
species was shown to have the same characteristics ob¬ 
served as a previously studied system investigating propi¬ 
onate |T] ; three steady-states and always stable. A second 
model introduced chlorophenol as the primary substrate 
and a dechlorinating organism forming a tri-culture that 
resulted in the emergence of bistable conditions between 
the complete washout steady-state and the two other vi¬ 
able states. Analysis of these bistabilities has shown that 
the desired operating condition (SS6) has the strongest 
basin of attraction such that complete washout can only 
occur when the initial condition for one of the biomass 
concentrations is zero. 

Although it is interesting to see the extent under 
which full chlorophenol mineralisation can occur using 
this standard model, of greater interest from an engi¬ 
neering perspective is the possibility of driving the sys¬ 
tem towards its limits of operational viability without 
compromising its function. With this in mind, the in¬ 
clusion of additional input terms for hydrogen and phe¬ 
nol was undertaken, and corresponding stability anal¬ 
ysis performed. In the case of hydrogen addition, it 
was expected that this would lead to a wider region 
of stability for the methanogen and, thus, extend the 
operational domain for SS6. Indeed, this was the case 
under relatively high concentrations of chlorophenol in¬ 
put, under which the methanogen population was main¬ 
tained up to a theoretical maximum, dictated by its max¬ 
imum growth rate. However at low chlorophenol input 
(< 0.5kgCOD/m 3 ), the nature of the resource competi¬ 
tion between the dechlorinator and methanogen is such 
that a switching of behaviour occurs in localities not con¬ 
ducive to stability of all three organisms. 

Under the standard three-tiered chlorophenol model, 
resource competition dictates that the dechlorinator is 
able to utilise the available hydrogen for growth more 
readily than the methanogen and under SS4 a compet¬ 
itive exclusion principle occurs and the latter is washed 
out. With hydrogen addition, however, an abundance 
of this resource allows the methanogen to access enough 


hydrogen to maintain a population at steady-state. How¬ 
ever, below a threshold concentration and under specific 
dilution rates, an excess of hydrogen leads to inhibition 
of the phenol degrader and washout. In effect, the ad¬ 
dition of hydrogen at this threshold results in a switch¬ 
ing between SS4 and SS5, rather than a stabilisation at 
SS6. The addition of phenol, however, under specific op¬ 
erational conditions can result in the emergence of SS6 
under low chlorophenol input conditions. Here, a spe¬ 
cific form of competition between the dechlorinator and 
methanogen occurs (as with SS6 in the standard model), 
by which both organisms benefit from the production 
of hydrogen by the intermediate phenol degrader, whilst 
the phenol degrader is stabilised by the syntrophic hy¬ 
drogen removal. In other words, the presence of the 
chlorophenol degrader allows for the production of phe¬ 
nol, which benefits both the phenol degrader and, indi¬ 
rectly, the methanogen. The presence of the methanogen 
reduces the inhibition on the phenol degrader, thus al¬ 
lowing the hydrogen resource to be maintained for both 
liydrogenotrophic populations. This can be seen as a 
form of mutualism rather than competition. The phenol 
addition plays a significant role under conditions where 
the chlorophenol degrader cannot produce enough phenol 
to maintain the phenol degraders, and under such condi¬ 
tions it is possible to achieve full chlorophenol mineralisa¬ 
tion beyond the standard model system. A schematic of 
these models and their interactions for the low chlorophe¬ 
nol input condition is shown in Fig. [8] 

Although higher dimensional systems are analytically 
restricted and lack generality, the work reported here for 
a three-tiered ’food web’ underlines the potential for ap¬ 
plying localised stability analysis within meaningful op¬ 
erating and parameter ranges, to identify properties of 
the system that both increase fundamental understand¬ 
ing of the behaviour of such microbial systems, but also 
guide thought on ways to manipulate or control them 
for potential process improvement. Further work will in¬ 
clude a thermodynamic rather than a kinetic inhibition 
term for the effect of hydrogen on the syntroph, and the 
extension of the model to polychlorinated phenols. 
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Appendix A: Parameter estimation 

The parameters of the model are based on the Chem¬ 
ical Oxygen Demand (COD) to be in consensus with 
the IWA Anaerobic Digestion Model No. 1 [5]. How¬ 
ever, only those related to hydrogen can be found in 
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*H2,ln 


using a conversion factor of 1 . 2 x 10 14 g/cell P 2 |, or 1 . 2 x 
10" 4 g/L. 

Therefore, the particulate CODx can be calculated by 


CODx= 2 gCOD/(g biomass)x 1.2x10 4 (g biomass)/L 
= 2.4x 10 -4 gCOD/L. 

The COD of the substrate per COD of biomass (fc m , c h) 
was then calculated using released chlor ions up to 
35/xM/Ld, which corresponds to COD consumed by 
biomass of 



FIG. 8. Schematic of observed microbial interactions in the 
three chlorophenol degradation models at low chlorophenol in¬ 
put concentrations. Top: Three-tier chlorophenol model with 
hydrogen inhibition of the phenol degrader, Middle: Hydro¬ 
gen addition model with phenol degrader washout. Bottom: 
Phenol addition model with extended stability. Organisms 
that are always present are shown as a green circle, a red 
circle is for when they are conditionally present, and a black 
circle represents biomass washout 


the ADM1 (for high-rate mesophilic processes). Subse¬ 
quently, the parameters for chlorophenol and phenol were 
chosen based on a combination of literature data and de¬ 
duction based on their chemical oxidation reaction. 

a. Parameters for chlorophenol Adrian et al. used 
the Dehalococcoides mccartyi strain CBDB1 to com¬ 
pletely convert 2,3-dichlorophenol, all six trichlorophe- 
nols, all three tetrachlorophenols, and pentachlorophenol 
to lower chlorinated phenols [TS], and their results are 
used partially here, with some deductive reasoning. The 
observed dechlorination rates in batch cultures with cell 
numbers of 10 7 mL _1 was found to be 35 /dVhD 1 . The 
COD of the CBDB1 strain was first estimated using the 
empirical formula for bacterial cell C 5 H 7 O 2 N, with cell 
numbers of 10'mL -1 . 

The oxidation reaction for C 5 H 7 O 2 N is 

C 5 H 7 O 2 N+ 5 O 2 -> 5C0 2 +NH 3 +2H 2 0. (Al) 

From this equation, it is easy to show that 1 mole 
of C 5 H 7 O 2 N requires 5 mole of O 2 , or 160gCOD. The 
molar mass of C 5 H 7 O 2 N is 113g/mol. Therefore, the 
proportion of C 5 H 7 O 2 N is 160/113 = 1.42g. By applying 
this to general organic matter, a factor of 2 is used. This 
means that each gram of biomass (e.g. strain CBDB1 
here) requires 2g of oxygen for growth. Subsequently, the 
cell number should be converted to grams (dry weight) 


35 /xmol[Cl ]/Ld 
2.4xl0 - 4 gCODx/L 


0.14mol[Cl ]/gCODxd. 


It should be noted that this result does not include 
monochlorophenol, which have to be deduced. For 
monochlorophenol, the COD must be converted from the 
COD of Cl - . By considering the stoichiometry of its oxy¬ 
genation, 


C 6 H 5 0Cl+6.50 2 -> 6C0 2 +2H 2 0+2H + +2Cr, (A 2 ) 

the COD of monochlorophenol can be calculated as 6.5x 
32 = 208gCOD/mol. From this, the maximum specific 
uptake rate for monochlorophenol degraders can be es¬ 
timated as k m ,ch = 0.14x208 = 29.12gCODs/gCODxd, 
where CODs is the substrate COD. 

For the half-saturation constant, Ks,ch, no data was 
found detailing empirical studies for monochlorophenol 
degradation. However, from a recent study on the anaer¬ 
obic microbial dechlorination of pentachlorophenol, some 
values were cited for dechlorinating bacteria using co¬ 
substrates of the 5-CP with other carbon sources [20]. A 
Ks value of 0.25/xM, which converts to 0.053kgCOD/m 3 
following the oxidation of monochlorophenol, was se¬ 
lected as a rough estimate for the chlorophenol degrader 
in this model. 

The yield constant, 1/h was calculated using infor¬ 
mation about the growth yield of the Dehalococcoides 
mccartyi strain CBDB1 on 2,3-dichlorophenol (DCP), in 
which 7.6xl0 13 cells per mole of CD were produced fT 8 j . 
Since yield is the ratio of CODx to CODs, it can be 
calculated by determining these components from their 
stoichiometry; 

: For 7.6xl0 13 cells, the dry weight is calculated 
as 7.6 x 10 13 x 1.2x lCD 14 g/cell = 0.912g. Know¬ 

ing that 2gCOD/gbiomass, stated earlier, then the 
COD x = 1.824g. 

CODs: One mole of CD released should be converted 
to COD that DCP consumed. It is known that: DCP—» 
2CD. That means one mole of CD comes from 0.5 mole 
of DCP. The oxidation equation of DCP is given by 
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C 6 H 4 0Cl2+60 2 -A 6C0 2 +H 2 0+2H++2Cr. (A 3 ) 

Therefore, 1 mole of DCP consumes 6 mole of oxygen, 
giving 6x32 = 192gCOD, which corresponds to 0.5 mole 
of DCP using 96gCOD. 

Subsequently, the yield of the chlorophenol degrader 
can now be calculated as 

Y ch = 1.824gCOD/96gCOD = 0.019gCOD x /gCOD s . 

The assumed half-saturation constant describing the 
affinity of the dechlorinator for hydrogen, I\s,h 2 ,c was 
arbitrarily set to lxlO -6 in the absence of supportable 
empirical data. The value was based on the known con¬ 
stant for the methanogen, and chosen to be an order of 
magnitude lower, accordingly. 

b. Parameters for phenol There is very limited data 
available on the kinetics of phenol degradation in anaero¬ 
bic environments. The parameters selected for the phenol 
degrader here are approximate values based on available 
literature and deduction from analogous substrate reac¬ 
tions. 

For example, the reaction of phenol is similar to that 
for propionate, but consumes twice as much oxygen. The 
oxidation equations for the two chemicals are written as 


stable. Therefore, by considering sets of parameters (D, 
•S'ch.in, 'S'ph.in and .SA^in) in turn on each steady-state, 
and through the use of Matlab, the complex polynomials 
for each steady-state in turn can be solved to determine 
which are meaningful and stable. Repeating this method 
with numerous sets of parameters leads to results show¬ 
ing where each steady-state is stable. 

In addition, some cases (SS3 or SS4) do not al¬ 
low straightforward expressions for the steady-state 
(quadratic and cubic equations for s o are obtained, re¬ 
spectively). In these situations, the above procedure 
is followed for all possible values of so at these steady- 
states. 


Appendix C: Stability analysis for the three-tier 
chlorophenol model with substrate addition 

Here, the expressions used for each steady-states in 
turn when hydrogen and phenol addition are included in 
the model are given. To obtain similar expressions when 
hydrogen and/or phenol addition is not present in the 
model, simply set Uh = 0 and/or u p = 0 , respectively. 
For the two-tier phenol model, consider sq = xq = 0 and 
solve accordingly. 

SSI (Xo = Xl = X2 = 0) 


Propionate : C 3 H 6 0 2 -|-3.50 2 —>■ 3C0 2 +3H 2 0 
Phenol : C 6 H 5 0H+70 2 -A 6C0 2 +3H 2 0. 

Hence, the activity of phenol is assumed to be dou¬ 
ble the activity of propionate. From the ADM1 [4], 
the maximum specific activity for propionate (/c m . pro ) 
is 13gCODs/gCODx, so k rn p h is assumed to be 
26gCOD s / gCOD x ■ 

The half-saturation constant for phenol was taken from 
a study by Eismann et al. 121 \ using a modified Hal¬ 
dane kinetic model from Edwards [22], which consid¬ 
ers the inhibitory effect of high substrate loading, for 
data fitting. At 35° C, the A'g, P h determined by the au¬ 
thors was 127mgPhenol/l, which converts to a value of 
0.302kgCOD/m 3 used in this model. 

In the absence of compelling literature data, the yield 
parameter for the phenol degraders, F p h, was taken to be 
equivalent to that for another acetogenic bacteria, the 
propionate degrader described in previous work [ 1 ] and 
having a value of 0.04gCODs/gCODx- 


Appendix B: Numerical Methods 


Equations (14), (15) and (16) result in so = uf, S\ = 
u g and s 2 = Uh, such that the steady-state is always 
meaningful. 

The Jacobian matrix, Jssi, can be written as 


—Ot+p, 0 — 

kA 0 

0 

0 

0 

0 ■ 


0 

—a+Mi ~kp 

0 

0 

0 

0 


0 

0 

-a+M 2 ~kc 

0 

0 

0 


-Mo 

0 

0 

— a 

0 

0 


woMo 

-Mi 

0 

0 

—a 

0 


—w 2 mo 

WiMi 

-M 2 

0 

0 

—a_ 


Its six eigenvalues are Ai 

= —a+Mo- 

kA, 

A 2 

= 

-a+ 

— kp, A 3 — 

—a+p2~kc, and A 4 = A 5 = 

- a 6 

= - 

-a. 

The 


condition for stability is Ai i2 ,3 < 0. By substituting so, 
si and s 2 into each of /xo, Mi and M 2 , stability is assured 
when 


u f 


Uh 


1 +Uf Kp+Uh 

4>1 u h 


—Ua < a; 


1 + Uh 


—kc < a. 


<j>l Ug 1 
1 + Ug 1 +KlUfi 


—kp < a; 


Although it is possible to simplify the system and ob¬ 
tain analytical expressions for characteristic polynomials 
for all of the steady-states solved in[C] it is impossible to 
gather the necessary information from these to determine 
the regions where each steady-state is meaningful and/or 


SS2 (xo = 0,Ki = 0, X 2 > 0) 


juations (14) and (15) give Sq = Uf and s 4 = u Q 
ing the remaining dimensionless equations 


i gives 
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a+kc 
a+kc-fo 
a ( u h -s 2 ) 
a+kc 


(Cl) 

(C2) 


Here, the Jacobian matrix, Jss 2 ? then becomes 


- a+no—kA 0 

0 — a+fii—k B 

0 0 

-Mo 0 

w oMo -Mi 

. -w 2 Mo wi/ii 

Its eigenvalues are Ai : 
Mi A 3 = A 4 = 


0 0 0 0 

0 0 0 0 

0 0 0 lx 2 

0 -a 0 0 

0 0 -a 0 

-(a+/cc) 0 0 —a—1x2. 

z —a+fi, q— fc^4, A2 = — a+ 
-Ix 2 +a 

-a, A 5j6 = -2- 1,2 

/2 

. For stability, it is sim¬ 


ply required that Ai 2 < 0 , as for all meaningful values of 
a , /, £2 and fc c , there is always As^^e < 0 . 


SS3: xo > 0, xi = 0, X 2 = 0 


Solving Eqs ( fl4| ), ( [TTj ) , ( [l5| and (16) gives 


SqW 2 (1 ~a-kA )+s 0 [w 2 Uf (a+A- 1) +Uh (1 - a -) 

— (a+A)( K Kp+u 2 )\ — (a+kA)(E P —oj 2 Uf+Uf l ) = 0 

(C3) 

S2 = U h -U 2 (Uf — So) (C4) 

a (u/-s 0 ) 


£0 = 


Si = M„+ 


a+fcA 
w 0 (o+Zca) 


(C5) 

(C 6 ) 


As mentioned previously in [B] when numerically check¬ 
ing whether this steady-state is meaningful, each possible 
value for so is considered in turn. 

The Jacobian matrix, Jss 3 > then becomes 


0 

0 

0 

0 

O 

Fx 0 

0 

-a — \- P i~k B 

0 

0 0 

0 

0 

0 

a+p 2 -kc 

0 0 

0 

~{a+k A ) 

0 

0 

1 

P 

1 

0 

O 

-Fx 0 

u 0 (a+k A ) 

-Mi 

0 

(jJqExq —a 

ui 0 Fx 0 

-u} 2 {a+k A ) 

Cdi/ii 

-M2 

—uj^Exq 0 - 

-a—u) 2 Fxo 

Its known eigenvalues 

are Ai 

= —G!+/il 

-k B , A 2 


—a+M 2 — kc, A 3 = —a. The other three eigenvalues are 
given by the characteristic polynomial 


(C7) 


where 


f 2 = 2a+Exo+Fu 2 x 0 (C 8 ) 

/1 = a 2 +kAExo+2Eaxo+kAFu}2Xo+2Fau} 2 xo (C9) 
/o = Ea 2 Xo+Fa 2 u> 2 Xo+kAEaXf)+kAFau} 2 Xo. (CIO) 

For stability it is r equired that Ai^ < 0 and all three 
roots from Eq. ( |C7[ ) have negative real parts. These are 
checked following the procedure described in [Bj 


SS4 (xo > 0, xi > 0, X 2 = 0) 


By solving Eqs. to ( fl(3| ) for £ 0 , x\, si and S 2 , and 
substituting these into equation (14), the following cubic 
expression for Sq is obtained 


P3S0+P2S0+P1S0+P0 = 0 (Cll) 


where 


P3 = ( kA+a) (( k B +a) (kA+a-I<iK P (kA+a)) 

-0i (k A +a )) {u h -WiUg+auf (w 2 -w 0 wi)) 

—K P (kA+a) (( k B +a) (kA+a—KjK P (kA+a)) 

-0i (kA+a))-u 1 (kA+a) ( k B +a) 

(kA+a-KjKp (k A +a)) (C12) 

P 2 = (((kB+a) (k A +a—KiK P (k A +a))-4>i (k A +a)) 
(kA+a —1)F(/l^iFq) ((k B +a) 

(kA+a—I<iK P (kA+a)— 1) 0i (kA+a— 1))) 
(uh-uiUg+auf (w 2 -w 0 wi))-A> (kA+a) 

((k B +a) (kA+a—KiK P (k A +a))-<fii (kA+a)) 

—K P (kA+a) ((k B +a) (k A +a—K I K P (kA+a)— 1) 

- 0i (kA+a-l))-a (kA+a) (w 2 -w 0 wi) ((k B +a) 
(k A +a-KiK P (k A +a))—(j)\ (k A +a)) 

-wi (kA+a) (k B +a) (kA+a—KjK P (kA+a)— 1) 

-wi (k B +a) (kA+a—1) (kA+a—I<iK P (kA+a)) 

(C13) 

Pi = ((k B +a) (kA+a—K I K P (fc^+a)-!) 

-0i (k A +a- 1)) (kA+a—1) ( u h -uiu g +au f (w 2 -w 0 w 1 )) 
—K P (kA+a) ((k B +a) (kA+a-KjK P (kA+a)— 1) 

-0i ( k A +a-l))-a (w 2 -w 0 wi) (((k B +a) 

(kA+a—I<iK P (kA+a)) 0i (kA+a)) (kA+a—1) 

-(- (kA+a) ((k B +a) (kA+a — Kjl\ p (kA+a) 1) 

- 0i (k A +a-l)))-Ui (k B +a) 

(kA+a— 1) (kA+a—KiK P (kA+a)— 1) (C14) 

Po = a (w 2 -w 0 wi) ((k B +a) (fc^+a-A'/A'p (fc^+a)-!) 
-0i (kA+a— 1)) (kA+a— 1). (C15) 


A 3 +/ 2 A 2 +/iA +/ 0 = 0 
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Once again, when numerically checking whether this 
steady-state is meaningful, each possible value for So 
in turn is checked. The remaining expressions for this 
steady-state are given by 


_ —Kp (a+k A ) (1+so) 
(l+so) (oi+k A )—so 
_ — (1 FKjS 2 ) {cx+kp) 

Sl (1+K i s 2 ) (a+k B )-(t> 1 
_ awi (u g -si)+a(u h -s 2 ) 
(a+k A ) (w 2 -woWi) 

_ uj 2 ( a+k A ) a 

Xl u}i(a+k B ) X ° LUi(a+k B ) 


(C16) 

(C17) 

(C18) 

(C19) 


The Jacobian matrix, Jgs 4 , becomes 


0 0 0 Ex 0 0 

0 0 0 0 Gx! 

0 0 — Q!+/X2 _ kc 0 0 

— (a+A:^) 0 0 —a—Ex 0 0 

uj 0 (a+kA) —(cn+fcs) 0 woExq —a—Gx 1 

-(jJ2{,OL+kA) uii(a+kB) — fJ-2 -UJ 2 EX 0 ujiGxi 


Fx 0 
Hx 1 
0 

-Fx 0 

UJqFXq — Hxi 
—ol—lo2Fxq+ojiHx\ 


= 3a+Exo+Gxi+Foj 2 xo—HujiXi (C21) 

fz = 3 a 2 +k A Exo+k B Gxi+ 3 Eaxo+ 3 Gaxi 
+k A Fu) 2 xo—k B HLJiXi+EGxoXi+3Fa(jj 2 Xo 
—3Hau>iXi+FGtj 2 XoXi—EHuJiXoXi 

—FGujqUiXqXi (C22) 

f 2 = a 3 +3Ea 2 x 0 +3Ga 2 Xi+3Fa 2 oj 2 x 0 —3Ha 2 uJiXi 
+2k A Eax 0 +2k B Gaxi+k A EGxoXi+k B EGxoXi 
+2k A Faoj 2 xo—2k B F[au}iXi+3EGaxoXi 
+k A FGui 2 XoXi—k A EF[uiiXoXi+k B FGu} 2 XoXi 
—k B EFluJiXoXi+3FGauj 2 xoXi—3EF[aijOiXoXi 
—3FGauiouJiXoXi—k A FGijJouJiXoXi 

—k B FGu)QijJiXQXi (C23) 

/1 = Ea 3 xo+Ga 3 xi+k A Ea 2 xo+k B Ga 2 xi+Fa 3 u] 2 xo 
-~Ha 3 bJiXi+k A Fa 2 Lo 2 xo—k B Ha 2 u)iXi+3kEGa 2 XoXi 
+3FGa 2 u} 2 xoXi—3EF[a 2 ujiXoXi+k A k B EGxoXi 
+2k A EGaxoXi+2k B EGaxoXi—3FGa 2 oJoUJiXoXi 
+k A k B FGu> 2 XoXi—k A k B EF[iOiXoXi+2k A FGaui 2 xoXi 
—2k A EF[au}iXoXi+2k B FGau} 2 XoXi—2k B EHaoJiXoXi 
—k A k B FGuioU)yxoXi—2k A FGauJoU]iXoXi 
—2k B FGau)QU)iXQXi (C24) 

/o = EGa 3 xoXi+FGa 3 ui 2 xoXi—EHa 3 ojiXoXi 
+k A EGa 2 XoXi+k B EGa 2 XoXi+k A FGa 2 co 2 XoXi 
—k A EF[a 2 uJiXoXi+k B FGa 2 Ljj 2 XoXi~k B EF[a 2 LOiXoXi 
—FGa 3 LOoU)iXoXi+k A k B EGaxoXi 
+k A k B FGaoj 2 XoXi—k A k B EF[au)iXoXi 
—k A FGa 2 u;oU)iXoXi—k B FGa 2 ujouJiXoXi 

—k A k B FGauiouJiXoXi. (C25) 


For stability it is required that Ai < 0 and all five roots 
from equation (C20) have negative real parts. 


SS5 (xo > 0,xi = 0,X2 > 0) 


Its known eigenvalue is Ai = —a+fi 2 —kc and the char¬ 
acteristic polynomial for the remaining five eigenvalues is 
given by 


By solving the dimensionless equations (14), (111, (151, 


(16) and (13), the steady-state is given by 


•^ 5 +/4^ 4 +/3^ 3 +/2^ 2 +/iA+/o — 0 (C20) 


_ -a-fcc 
a+kc~4> 2 

_ — (Kp+s 2 ) ( a+k A ) 
(Kp+s 2 ) {cx+k A )—s 2 
a(u f -s 0 ) 

a. (u h -s 2 )-u 2 ( a+k A ) x 0 
X2 ~ a+k c 

ix 0 ( a+k A )x 0 

Sl = Ug ~|-. 

a 


(C26) 

(C27) 

(C28) 

(C29) 

(C30) 


where 
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The Jacobian matrix, Jss 5 ; becomes 


0 

0 

0 

Ex o 

0 

Fx 0 

0 

-a+ni—ks 

0 

0 

0 

0 

0 

0 

0 

0 

0 

lx 2 

a+kA 

0 

0 

—a—Ex 0 

0 

-Fx 0 

ui 0 ( a+k A ) 

-Mi 

0 

Eui 0 x 0 

—a 

FuIqXq 

_-w 2 ( a+k A ) 

wiMi 

- (a+kc) 

—Ew 2 x 0 

0 

-a—Fu2Xo—Ix2_ 


Its known eigenvalues are Ai = — a+ni —k B and A 2 = 
— a. The remaining four are given by the characteristic 
polynomial 


A 4 +/3A 3 +/ 2 A 2 +/!A+/o = 0 (C31) 

where 


y*3 — 2cx-\-Exq-\-Ix2~\~Fuj2Xq (C32) 

f-2 = a 2 +kAEx 0 +CIx2+2Eax 0 +2IaX2+kAFuj2Xo+ 
EIxoX2+2Fauj2Xo (C33) 

fi = Ea 2 xo+Ia 2 X2+Fa 2 uj2Xo+kAEaxo+kcIax2 
+kAEIxoX2+kcEIxoX2+kAFaco2Xo+2EIaxoX2 

(C34) 

/o = EIa 2 XoX2+kACEIxQX2+kAEIaxoX2+kcEIaxoX2- 

(C35) 


For stability it is req uired that Ap 2 < 0 and all four roots 
from equation (C31|) have negative real parts. 


SS6 (xo > 0 , 2:1 > 0, X 2 > 0) 


Solving all six dimensionless equations to ( |T6| ), 
leads to the expressions for this steady-state 


s 2 

si 

so 

*0 

Xi 

X2 


-a-kc 

Ot+kc~4>2 

— (1+A/s 2 ) (a+k B ) 
(1+A'/s 2 ) (a+k B )-(j> 1 

- (K P +s 2 ) ( a+k A ) 
(Kp+s 2 ) {a+k A )-<t> 1 
a (uf-sp) 

(a+k A ) 

a (Mg-si)+wo (a+kA) %o 
a+k B 


(C36) 

(C37) 

(C38) 

(C39) 

(C40) 


a (uh—s 2 )—W 2 (a+kA) a-’o+wi (a+k B ) X\ 
a+kc 


(C41) 


The Jacobian matrix, Jss6, becomes 


0 

0 

0 

Ex 0 

0 

Fx 0 

0 

0 

0 

0 

Gx 1 

Hx 1 

0 

0 

0 

0 

0 

1X2 

~(a+k A ) 

0 

0 

O 

H 

1 

d 

1 

0 

-Fx 0 

u 0 (a+k A ) 

~(a+k B ) 

0 

ui 0 Ex 0 

—a—Gx 1 

UJqFXq — HXi 

-u 2 (a+k A ) u 1 (a+k B ) 

-(a+kc) 

-ui 2 Ex q 

uj\Gx\ 

—a.—uj 2 Fxo+uJiF[xi 


Here, all six eigenvalues are given by the characteristic 
polynomial 

^ 6 +/5^ 5 +/4A 4 +/3A 3 +/ 2 A 2 +/iA-|-/o = 0 (C42) 

where 

/5 = 3a+Exo+Gxi+Ix2+Foj2Xo—HujiXi (C43) 

fi - 3a 2 +kAFxQ+k B Gx\+k(jIx2+3EaxQ+3Gax\ 
+3Iax2+kAFaj2Xo—k B HujiXi+EGxoXi+EIxoX2 
+GIx\X2+3Fau>2XQ—3HauJiXi+FGuj2XQXi 
-~EHliJiXqX\—FGujqujiXoXi (C44) 

/j = a 3 +3Ea 2 xo+3Ga 2 Xi+3Ia 2 X2+3Fa 2 uj2Xo 
—3Ha 2 uiiXi+2kAEaxo+2k B Gaxi+2kcIaX2 
+kAEGxoXi+k B EGxoXi+kAEIxoX2+kcEIxoX2 
+k B GIxiX2+kcGIxiX2+2kAFau}2Xo—2k B Hau!iXi 
+3EGaxoXi+3EIaxoX2+3GIaxiX2+kAFGtU2XoXi 
~kAEHu}iXoXi+k B FGu}2XoXi—k B EHujiXoXi 
+EGIxoXiX2+3FGau>2XoXi—3EHau>iXoXi 
—3FGauiouJiXoXi—kAFGtjJoU)iXoXi 

—ksFGcooUJiXoXi (C45) 

/ 2 = Ea 3 xo+Ga 3 xi+Ia 3 X2+kAEa 2 xo+k B Ga 2 xi 
+kc la 2 X 2 +Fa 3 1 x 2 X 0 —Ha 3 LOiXi+kAFa 2 1 x 2 X 0 
—k B Ha 2 uJiXi+3EGa 2 xoXi+3EIa 2 XoX2+3GIa 2 XiX2 
+3FGa 2 iX2XoX\—3EF[a 2 LX\XoXi+kAk B EGxoXi 
+kAkcEIxoX2+k B kcGIx\X2+2kAEGaxoXi 
+2k B EGaxoX\+2kAEIaxoX2+2kcEIaxoX2 
+2k B GIaxiX2+2kcGIaxiX2—3FGa 2 ixo'XiXoXi 
+kAk B FGu>2XoXi~kAk B EF[LOxXoXi+kAEGIxoXiX2 
+k B EGIxoXiX2+kcEGIxoXiX2+2kAFGauj2XoXi 
—2kAEF[au>iXoXi+2k B FGau>2XoX\—2k B EF[au>iXoXi 
+3EGIaxoXiX2—kAk B FGujoixiXoXi—2kAFGaujoLOiXoXi 

—2k B FGaojouJiXoXi (C46) 

/1 = Ea 3 xo+Ga 3 Xi+kAEa 2 xo+k B Ga 2 xi+Fa 3 1 x 2 X 0 

—Ha 3 u)iXi+kAFa 2 L02Xo—k B Ha 2 u)iXi+3EGa 2 xoXi 

+3FGa 2 oj2XoXi—3EF[a 2 ujiXoXi+kAk B EGxoXi 

+2kAEGaxoXi+2k B EGaxoXi—3FGa 2 ixoUJiXoXi 

+kAk B FGu>2XoXi~kAk B EF[iiJxXoXi+2kAFGaijJ2XoXi 

—2kAEF[auiiXoXi+2k B FGaui2XoX\—2k B EHauiiXoXi 

-kAk B FGujoLOiXoXi—2kAFGaujoUJ-LXoXi 

—2k B FGaujou>iXoXi (C47) 

fo = EGIa 3 XoXiX 2 +kAEGIa 2 xoX\X 2 +k B EGIa 2 xoX\X 2 

+kcEGIa 2 xoX\X2+kAk B kcEGIxoXiX2 
+kAk B EGIaxoX\X2+kAkcEGIaxoX\X2 
+k B kcEGIaxoX\X2■ (C48) 
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For stability it is required that all six roots from equa¬ 
tion (C42) have negative real parts. 


SS7 (xo = 0 , 2:1 > 0 , 2:2 = 0) 


Solving equations (141, (15), (12) and (16) gives 


s\ [KjUh ( a+k B )-a-k B +K I u!iUh ( a+k B )] + 
si [Kiu h (a+k B )-a-(t>i-Kiu}i ( a+k B )-k B 
+ KiuiUh {a+k B )]—KiUi ( a+k B ) = 0 

__ 4>isi _i_ 

S2 (1 +s 1 )K I (a+k B ) /\ 7 

S 0 = Uf 

a ( Ug-si ) 

Xl ~ a+k B ' 


(C49) 

(C50) 

(C51) 

(C52) 


As mentioned previously in[Bj when numerically check¬ 
ing whether this steady-state is meaningful, each possible 
value for s 0 is considered in turn. 

The Jacobian matrix, Jss 7 , becomes 


—a—kc 
cx+kc~4> 2 

_ — (I+A/S 2 ) ( a-\-k B ) 

Sl ( 1 +K iS2 ) (a+k B )—<f>i 
a (ug-sx) 

Xl = ^kT 

a (%-s 2 )+wi ( a+k B ) x\ 
X2 ~ a-+k t c 

So = Uf 


(C57) 

(C58) 

(C59) 

(C60) 

(C61) 


The Jacobian matrix, Jsssj becomes 


—a+Ho—kA 

0 

0 

0 

0 

0 

0 

0 

0 

0 

GX! 

Hx 1 

0 

0 

0 

0 

0 

1 X 2 


0 

0 

—a 

0 

0 

CJ 0 M 0 

- ( a+k B ) 

0 

0 

—a—Gx 1 

—Hx 1 

—^2 Mo 

wi ( a+k B ) - 

(a+kc) 

0 

Glo\X\ 

— a.+Hu)\X\ —lx 2_ 


O+Mo — kA 

0 

0 

0 

0 

0 

0 

0 

0 

0 

Gx 1 

Hx 1 

0 

0 

-a+^ 2 -k c 

0 

0 

0 

Mo 

0 

0 

—a 

0 

0 

UlofJ .0 

- (a+fcs) 

0 

0 

—a—Gx 1 

— Hx 1 

— oi 2 Mo 

wi (a+A: s ) 

-M 2 

0 

Gui\ x 1 

—a+Hui 1X1. 


Its known eigenvalues are Ai = — a+no~kA and A 2 = 
—a. The remaining four eigenvalues are given by the 
characteristic polynomial 


A 4 +/3A 3 +/2A 2 +/!A+/o = 0 (C62) 


Its known eigenvalues are Ai = —a+no—hA, A 2 -I— a+ 
M 2 ~kc and A 3 = —a. The other three eigenvalues are 
given by the characteristic polynomial 


where 


•^ 3 +/ 2 A 2 +/iA-|-/o — 0 (C53) 

where 


f 2 = 2a+Gx\—HujiX\ (C54) 

fi = a 2 +k B Gxi+2Gaxx—k B HujiXi—2HcdjJiXi (C55) 
fo = Ga 2 xi—Ha 2 uj\X\+k B Gax\—k B HaujiXi. (C56) 


For stability we re quire Ai ;2 ,3 < 0 and all three roots 
from equation (C531 to have negative real parts. These 
are checked following the procedure described in [B| 


/3 = 2a+Gx\+Ix 2 —HuiiX\ (C63) 

f 2 = a 2 +k B Gxi+kc!x 2 +2Gaxi+2Iax 2 —k B HujiXi+ 
GIx\x 2 —2Ha<jJiXi (C64) 

fi = Ga 2 Xi+Ia 2 x 2 —Ha 2 coiXi+k B Gaxi+kcIcxx 2 
+k B GIxix 2 +kcGIxix 2 ~k B HaujiXi 
+2GIax\X 2 (C65) 

fo = GIa 2 xix 2 +k B kcGIxix 2 +k B GIaxix 2 
+kcGIax\x 2 . (C 66 ) 


For stability it is required that Ai .2 < 0 and all four 
roots from equation (C62) to have negative real parts. 


SS8 ( 2:0 = 0 , 2:1 > 0 , 2:2 > 0) 


By solving the dimensionless equations (141, (151, (12), 


(16) and (13), the steady-state is given by 


REFERENCES 
























19 


[1] A. Xu, J. Dolfing, T. P. Curtis, G. Montague, E. Martin, 
Maintenance affects the stability of a two-tiered microbial 
‘food chain’?, J. Theor. Biol. 276 (1) (2011) 35-41. doi: 
10.1016/j.jtbi.2011.01.026 

[2] J. Dolfing, J. M. Tiedje, Hydrogen cycling in a three¬ 
tiered food web growing on the methanogenic conversion 
of 3-chlorobenzoate, FEMS Microbiol. Ecol. 38 (5) (1986) 
293-298. doi : 10.1111/j.1574-6968.1986.tb01740.x 

[3] J. Dolling, J. M. Tiedje, Kinetics of two complementary 
hydrogen sink reactions in a defined 3-chlorobenzoate de¬ 
grading methanogenic co-culture, FEMS Microbiol. Ecol. 
86 (1) (1991) 25-32. doi:10.1111/j.1574-6968.1991. 
tb04792.x 

[4] D. J. Batstone, J. Keller, I. Angelidaki, S. V. Kalyuzh¬ 
nyi, S. G. Pavlostathis, A. Rozzi, W. T. M. Sanders, 
H. Siegrist, V. A. Vavilin, Anaerobic Digestion Model 
No. 1, Tech. Rep. Report No. 13, IWA Publishing, Lon¬ 
don, UK (2002). 

[5] O. Bernard, M. Polit, Z. Hadj-Sadok, M. Pengov, 
D. Dochain, M. Estaben, P. Labat, Advanced monitor¬ 
ing and control of anaerobic wastewater treatment plants: 
software sensors and controllers for an anaerobic digester, 
Wat. Sci. Technol. 43 (7) (2001) 175-182. 

[6] D. Gaida, C. Wolf, C. Meyer, A. Stuhlsatz, J. Uppel, 
T. Back, M. Bongards, S. McLoone, State estimation for 
anaerobic digesters using the adml, Wat. Sci. Technol. 
66 (5) (2012) 1088-1094. doi:10.2166/wst.2012.286 

[7] V. A. Vavilin, S. V. Rytov, L. Y. Lokshina, J. A. Rintala, 
G. Lyberatos, Simplified hydrolysis models for the opti¬ 
mal design of two-stage anaerobic digestion, Wat. Res. 
35 (2001) 4247-4251. doi:10.1016/s0043-1354(01) 
00148-8 

[8] B. Benyahia, T. Sari, B. Cherki, J. Harmand, Bifurcation 
and stability analysis of a two step model for monitor¬ 
ing anaerobic digestion processes, J. Proc. Control 22 (6) 
(2012) 1008-1019. doi:10.1016/j.jprocont.2012.04. 
012 

[9] S. Shen, G. C. Premier, A. Guwy, R. Dinsdale, Bi¬ 
furcation and stability analysis of an anaerobic diges¬ 
tion model, Nonlinear Dynam. 48 (4) (2007) 391-408. 
doi:10.1007/s11071-006-9093-1 

[10] I. Simeonov, S. Diop, Stability analysis of some nonlinear 
anaerobic digestion models, Int. J. Bioautomation 14 (1) 
(2010) 37-48. 

[11] B. C. Baltzis, A. G. Frederickson, Coexistence of two 
microbial populations competing for a renewable resource 
in a non-predator-prey system, Bulletin of Math. Biol. 46 


(1984) 155-174. doi:10.1007/bf02463728 

[12] R. M. May, Stability and Complexity in Model Ecosys¬ 
tems, Princeton University Press, 1973. 

[131 S. L. Pimm, Food Webs, University of Chicago Press, 
Chicago, 2002. 

[14] B. Dubey, J. Hussain, Modelling the interaction of two 
biological species in a polluted environment, J. Math. 
Anal. Appl. 246 (1) (2000) 58-79. doi: 10.1006/jmaa. 
2000.6741 

[15] J. Wu, H. Nie, G. S. K. Wolkowicz, A mathematical 
model of competition for two essential resources in the 
unstirred chemostat, Siam J. Appl. Math. 65 (1) (2004) 
209-229. doi : 10.1137/s0036139903423285 

[16] A. Bornhoft, R. Hanke-Rauschenbach, K. Sundmacher, 
Steady-state analysis of the Anaerobic Digestion Model 
No. 1 (ADM1), Nonlinear Dynam. 73 (1-2) (2013) 535- 
549. doi:10.1007/sll0710130807x 

[17] M. Weedermann, G. Seo, G. S. K. Wolkowicz, Mathemat¬ 
ical model of anaerobic digestion in a chemostat: Effects 
of syntrophy and inhibition, J. Biol. Dynam. 7 (1) (2013) 
59-85. doi:10.1080/17513758.2012.755573 

[18] L. H. Adrian, S. K. Hansen, J. M. Fung, H. Gorisch, 
S. H. Zinder, Growth of delialococcoides strains with 
chlorophenols as electron acceptors, Env. Sci. Technol. 
41 (7) (2007) 2318-2323. doi:10.1021/es062076m 

[19] F. E. Loftier, J. Yan, K. M. Ritalahti, L. Adrian, E. A. 
Edwards, K. T. Konstantinidis, J. A. Muller, H. Fuller¬ 
ton, S. H. Zinder, A. M. Spormann, Dehalococcoides 
mccartyi gen. nov., sp. nov. obligately organohalide- 
respiring anaerobic bacteria relevant to halogen cycling 
and bioremediation, belong to a novel bacterial class, De- 
halococcoidia classis nov., order Dehalococcoidales ord. 
nov. and family Dehalococcoidaceae fam. nov., within the 
phylum chloroflexi, Int. J. Syst. Evol. Microbiol. 63 (2) 
(2012) 625-635. doi:10.1099/ijs.0.034926-0 

[20] Z. Li, Y. Inoue, S. Yang, N. Yoshida, A. Katayama, 
Mass balance and kinetic analysis of anaerobic micro¬ 
bial dechlorination of pentachlorophenol in a continuous 
flow column, J. Biosci. Bioeng. 110 (3) (2010) 326-332. 
doi: 10.1016/j.jbiosc.2010.03.008 

[21] F. Eismann, P. Kuschek, U. Stottmeister, Microbial 
phenol degradation of organic compounds in natural 
systems: Temperature-inhibition relationships, Environ. 
Sci. Pollut. Res. Int. 4 (4) (1997) 203-207. doi: 10.1007/ 
bf02986346 

[22] V. H. Edwards, The influence of high substrate concen¬ 
trations on microbial kinetics, Biotechnol. Bioeng. 12 (5) 
(1970) 679-712. doi:10.1002/bit.260120504 





